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ABSTRACT 

We place constraints on the formation redshifts for blue globular clusters 

(BGCs), independent of the details of hydrodynamics and population III star 
formation. The observed radial distribution of BGCs in the Milky Way Galaxy 
suggests that they formed in biased dark matter halos at high rcdshift. As a 
result, simulations of a ~ 1 Mpc box up to 2; ~ 10 must resolve BGC formation 
in ACDM. We find that most halo stars could be produced from destroyed BGCs 
and other low- mass clusters that formed at high redshift. We present a proof-of- 
concept simulation that captures the formation of globular-like star clusters. 

Subject headings: early universe - globular clusters: general - galaxies: star 
clusters 



1. Introduction 

The first stars and star clusters are responsible for the initial chemical enrichment of 
the universe. They are signatures of the earliest stages of galaxy formation, and their rem- 
nants may provide the seeds for supermassive black holes. The oldest known star clusters 
arc globular clusters (GCs). In the Milky Way, their color distribution is bimodal (Zinn 
1985; Harris 2009), possibly indicating two formation channels. The blue, metal-poor GCs 
([Fe/H] < — 1) have a typical [Fe/H] ~ — 1.5 and the red, metal-richer ([Fe/H] > — 1]) pop- 
ulation has a typical [Fe/H] ~ —0.5. Blue GCs (BGCs) have chemical signatures and radial 
distributions that are similar to halo stars (e.g., Helmi 2008), while the red population is 
associated with the galactic disk and shows a clear circular velocity component. Relative 
age estimates show that red GCs are younger than BGCs, with an average age separation 
of about 1.5 Gyr between the Galactic populations (e.g., de Angeh et al. 2005). The typical 
GC in the Galaxy has an My ~ —7.3 (Harris 1991) and mass M — 1.4: x lO^M©, assuming 
a mass-to-light ratio Ty = 2. The dark matter content of these systems is small, if present 
at all (e.g., Moore 1996; Baumgardt et al. 2009). 
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The formation of GCs remains poorly understood (see Brodie & Strader [2006] for a 
summary). Motivated by their extreme ages, Peebles & Dicke (1968) argued for their forma- 
tion shortly after recombination, preceding protogalaxy formation. The GCs are centrally 
concentrated in the galaxy with half of them closer to the galactic center than the Sun. In 
contrast, the half mass radius of the dark halo is over 10 times larger. As a result, most 
subsequent theoretical work has focused on producing them after the galaxy collapses. These 
ideas have included compression behind strong shocks (Gunn 1980) or thermal instabilities 
is the early protogalaxy (Fall & Rees 1985). However, these conditions could never have 
existed in dwarf galaxies, yet some contain BGCs (e.g., Fornax and WLM; Harris 1991). So 
while we reject these ideas here, they can be revived to produce the younger population of 
GCs in mergers (Ashman & Zepf 2001). 

BGC ages remain a compelling case for formation before the collapse of the galaxy 
(e.g., Sarajedini et al. 2007). In this case, an alternative to explain the radial concentration 
is biasing owing to their early formation in rare peaks. Moore et al. (2006) showed that 
the radial distribution of dark matter that collapsed into halos by 2; ~ 12, with mass scales 
> 2.5a{M, z), is sufficiently biased to match the radial distribution of the Galaxy's BGCs 
and halo stars. A few simulations have tried to capture this early formation in dwarf galaxies. 
Bromm & Clarke (2002) identified gas clumps with masses ~ lO^M© in their simulation of 
a proto-dwarf galaxies (~ 10^ Mq) as possible proto-GCs, while Kravtsov & Gnedin (2005) 
used a a subgrid model to split "cluster particles" in dwarf galaxies with M > IO^Mq. 
However, these simulations have not yet demonstrated that GCs do form, nor why a lO^Mo 
cloud should form a GC at high z, while such clouds form associations or open clusters today. 

Halo stars are also old, but their connection with BGC formation remains uncertain. 
These stars have abundances similar to BGCs and a consistent radial distribution (e.g., 
Brodie & Strader 2006; Helmi 2008). In the Galaxy, the cluster mortality rate is high (Lada 
& Lada 2003; Bastian & Goodwin 2006), and the destruction of BGC-like star clusters that 
formed could produce a substantial fraction, even all, of the halo. Fall & Zhang (2001) 
showed that a power law initial cluster luminosity function (ICLF), like what is seen in 
the Antennae (Schweizer 1987; Zhang & Fall 1999), could be shaped into the present-day 
globular cluster luminosity function (GCLF). 

We focus on BGC formation using arguments that are independent of the details of 
hydrodynamics and star formation. We show that halo stars and BGCs could have a common 
origin; we place constraints on the epoch of halo star and BGC formation; and, we argue 
that simulations of a biased region of the universe, with box sizes ~ 1 Mpc comoving, should 
capture BGC formation by 2; ~ 10, regardless of the mechanism. In §2, we estimate the mass 
of the halo star population in the Milky Way Galaxy (MWG) by combining the GCLF with 
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an ICLF. We then find the maximum redshift that BGC and halo star formation could have 
been completed in §3. In §4, wc use radial biasing to constrain star formation efficiencies. 
Wc present in §5 a proof-of-concept simulation of the formation of star clusters at high z. 
Our conclusions are summarized in §6. 

2. Luminosity Functions and Mass Estimates 

In this section, we calculate the total initial mass of MWG halo stars and BGCs if 
they formed with an Antennae-like initial cluster mass function (ICMF). The GCLF is 
$ oc exp(— — Vo)^/2(Ty), for absolute visual magnitude V and turnover magnitude Vq. The 
globular cluster mass function (GCMF) can be estimated by assuming a constant mass-to- 
light ratio Ty: /qc oc cxp(— (yU — yUo)^/2(T^), where ^ = logm for mass m. Wc set Vq = —7.3 
and ay = 1.2, which is consistent with the Harris (1991) and Gnedin (1997) fits for MWG 
GCs. Assuming Ty = 2 gives Hq = 5.15 and a^ — OAS. We normalize the mass function to 
the Milky Way's current BGC population adopting a total number Nqc ~ 100. To represent 
the ICMF, we follow Zhang & Fall (1999) with an ICMF /5c = mdN/dm oc m-\ Figure 

1 shows curves for fee and fsc, where fsc is normalized such that it osculates fee at 
A*osc = /^o + 2.3(j^. Integrating fsc yields a total mass of 104mosc (A*high — /^low)- Including 
star clusters between -2 and -12 V, the total mass in the initial power law distribution is 

2 X 10^ Mq, compared with ~ 2 x 10'^ M© for the present-day BGCs. 

This calculation underestimates the ICMF as even clusters with masses of mosc will 
suffer destruction. Typical GC destruction rates are 0.03/Gyr arc 0.1/Gyr for evaporation 
and evaporation+disk-shocking, respectively (Gnedin & Ostriker 1997). Assuming these are 
constant over 13 Gyr, the mass function is reduced by a factor of 1.5 or 4. Therefore the mass 
in halo stars that could be produced during BGC formation is between 3 and 8 x lO^M©, 
consistent with estimates of the old stellar halo mass (e.g., Binney & Tremaine 2008). 

Other authors have investigated whether the halo population could be produced by 
dissolved clusters. Surdin (1995) concluded that cluster dissolution could not produce the 
mass the entire spheroid (bulge+halo; ~ 5 x lO^M©) (see also Ostriker & Gnedin 1997). 
However, the bulge population is chemically different from the stellar halo and the BGCs. 
The peak of the metaUicity distribution in the bulge is [Fe/H] ~ —0.1 and extends to ~ — 1 
(Zoccah et al. 2003), while the halo and BGCs have peaks at [Fe/H] ~ —1.6 and extend to 
~ —3 (Helmi 2008, Zinn 1985). Present thoughts on the bulge are that it is a pseudo-bulge 
formed through disk processes (see Binney 2009 for a review) and should be left out of this 
comparison. 
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Henceforth, we assume that the halo star population came from dissolved clusters that 
formed along with the present-day BGCs. It is unclear whether the Antennae-like ICMF 
describes star clusters at high z, but the detailed shape of the ICMF is not central to 
our argument. We only require that a distribution of clusters is produced along with the 
BGCs, and an Antennae ICMF gives us a fiducial model. We find a minimal extrapolated 
mass of 2 X lO^M©. Including cluster evaporation, the mass in halo stars produced during 
BGC formation is consistent with 3-8 x IO^Mq. Hereafter, we assume this mass to be 
Mh ~ 5 X lO^M©, consistent with observational constraints on the MWG's stellar halo. 



3. Required Star Formation Efficiency 

We now estimate the maximum rcdshift that BGC formation could have been completed 
by requiring enough baryons in appropriate collapsed objects to form the old stellar halo. 
To do so, we use the Sheth & Tormen (2002) formahsm for the conditional mass function. 
We set the critical density for collapse Sc — 1.686. The value for a density fluctuation 5 
for a given z is determined by scaling 5c with the growth function. The mass variance is 
calculated using a top-hat fllter in real space. We use the Eisenstein & Hu (1998) transfer 
function with baryonic acoustic oscillations. The conditional mass function is integrated to 
find the mass in halos with M > M m\-n at 2; in a halo of mass Mq at Zq. 

We take Mq = IO^^Mq at zq = 0, and assume WMAP5 cosmology (Qc — 0.214, 
I). I) I ;. = 0.742, and h = 0.719; Hinshaw et al. 2009). The slope of the power spectrum 
n — 0.95, and erg = 0.8. The minimum mass for a halo M„^i^ is set by the lowest mass we 
expect to be relevant to BGC formation. This mass must be > 10^ Mq, i.e., where we expect 
the first stars to form (e.g., Bromm et al. 2009). Figure 2 shows the cumulative baryon mass 
in collapsed objects, assuming the universal baryon fraction, for minimum dark matter halo 
masses of 10^, 10*^, 10*, and IO^Mq. The two thin horizontal fines show the mass fraction 
in baryons required to form ~ 5 x IO^Mq of halo stars and BGCs, assuming a global star 
formation efficiency (SFE) of 77 = 1% (top) and 10% (bottom). Reionization may well be 
responsible for the truncation of BGC formation (see Brodie Sz Strader 2006), hence we show 
the WMAP5 reionization redshift with gray shading in Figure 2. The plot shows that halo 
star and BGC formation could have been completed before reionization for SFEs ~ 10%, if 
BGC formation took place in halos less massive than IO^A/q. This does not mean that BGCs 
could not form in higher mass halos; instead, it means that they could not form exclusively 
in high-mass halos. Although reionization could be responsible for halting BGC formation 
after 2; ~ 10, reionization is not central to our argument, as we shall see next. 
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4. Radial Biasing Constraints 

Diemand et al. (2005) and Moore et al. (2006) used dark matter only simulations to 
demonstrate that material that collapses into halos with mass scales of z/ > 2.5 at high 
redshift, where v = 5c/(j{M, z), has a radial distribution in the z = galaxy that is consistent 
with the MWG halo stars' and BGCs'. We note that the characteristic mass, i.e., the knee in 
the Press-Schechter function, corresponds to u — ^/6 (Binney & Tremaine 2008). Using the 
Sheth-Tormen conditional mass function, we integrate over all halos with masses greater than 
(t(M, z) for 1/ = 2.5. At ^ ~ 12, this mass fraction fu>2.5 ~ 0.027 (thick solid line in Figure 
2). To produce the halo and BGC population from exclusively high-z/ halos, the global SFE 
needs to be > 10%. We analyzed the Diemand et al. simulation data for comparison with 
the Sheth-Tormen formalism. We reproduce the biasing results, so the radial distributions of 
halo stars and BGCs are consistent with u = 2.5 halos from z > 10. However, more collapsed 
objects were found in the simulations than predicted by the conditional mass function for 
z >5. The biased mass calculated from Sheth-Tormen may be a lower hmit. 

For a 10% SFE, about 10% of halo stars and BGCs could be produced by z ~ 18. For the 
same SFE, the rest of the stars can be produced by 2; ~ 13, 120 Myr later. The Lagrangian 
volume for the MWG is approximately 27 Mpc^ comoving. Using these figures, we estimate 
the star formation rate required to form 5 x IO^Mq of stars is ~ 0.15Mq yr~^ Mpc~^. This 
value is similar to the observed peak in the universal star formation rate (Madau et al. 1998; 
0.12-0.17 Mq yr~^ Mpc~^ at z ~ 1.5). If BGC formation extends to z ~ 10, our estimate 
would be about a factor of two lower. Note that the universe is only about 500 Myr old by 
z ~ 10; de Angeli et al. (2005) found an age spread in their low-metallicity GC sample that 
is less than 600 Myr and consistent with zero. 

Combining our above arguments, we conclude that capturing the formation of halo 
stars and BGCs in ACDM using numerical techniques requires only a small, biased box. 
BGCs should form in a ~ 1 Mpc comoving box by 2; ~ 10, regardless of the underlying 
process. Based on the Sheth-Tormen mass function, the simulation volume should contain 
three lO^M© halos by the end redshift. If BGCs do not form in such a simulation, there is 
likely a problem with the subgrid model or the underlying cosmology. 
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5. Test Simulation 
5.1. Setup 

We highlight preliminary simulation results that will be presented in detail in a forth- 
coming paper. We employ RAMSES (Teyssier 2002), an adaptive mesh refinement cosmol- 
ogy code that solves the equations of hydrodynamics using a second-order unsplit Godunov 
method and integrates particles using particle-mesh techniques. The code has been aug- 
mented to follow nonequilibrium chemistry for a nine-species gas: e, HI, Hll, Hel, Hell, 
Helll, H^, H2, and H^. Cooling due to collisions between H2 and the other species is mod- 
eled according to Glover & Abel (2008). The Dubois & Teyssier (2008) star formation recipe 
is modified to capture pop 111 star formation and supernova feedback. Wherever the gas 
number density is greater than the critical value Uq — lO^cm""^, there is a chance to form a 
star particle. The local SFE is set to 1% per dynamical time. 

If the metallicity of the gas is above 2'thresh = 5 x 10~'^Zq, the resulting star particle 
represents a distribution of population II stars. This threshold is lower than suggested by 
Bromm & Loeb (2003) because we wanted to explore whether very low-metallicity population 
II stars could form if given the opportunity. Ten percent of the star particle's mass is injected 
back into medium after 10 Myr as supernova ejecta, with 0.1 of the ejecta mass in metals. 
Below ^thresh, the particles represent a distribution of population III stars. 

The IMF for population III stars remains unclear. The minimum mass may be > IOMq 
(e.g, Bromm et al. 2009), but this remains highly uncertain. Likewise, we know neither 
the maximum mass for a pop III star nor the importance of pair-instability supernovae 
to the enrichment of the second generation of stars. Tumlinson et al. (2004) argue that 
the observed element ratios of old pop 11 stars are inconsistent with enrichment by pair- 
instability supernovae. Moreover, the IMF shows little sensitivity to metallicity (Kroupa 
2001). Whether this extends to extremely low metallicities remains unknown. In the absence 
of better constraints on pop 111 star formation, we assume a Salpeter IMF (dN/dM ~ 
^-2.35'^ with a low-mass cutoff M^o (see below). Whenever a star particle is formed, the 
particle's mass Mgp is distributed among a random distribution of stars. The probability 
that a star is between masses Mq and Mi is given by P = A(Mf^'^^ — M^^'^^), where 
A = (Mtop^^ - M"^-^^)-^ and Mtop = IOOOMq is the most massive star in the distribution. 
A given star mass is selected by setting Mgtar = {r/A + M~^-^^)~^/^-^^, where r is a random 
variable with an even distribution between and 1. This process continues until the total 
mass in stars is equal to Mgp. We set M^o = 3Mq, which is based on the theoretical 
calculations of Nakamura & Umemura (2001), who find that the smallest Jeans mass for 
which the gas becomes opaque to H2 lines is ~ IM©. 
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For each population III star in the distribution, its contribution to ejecta, metal enrich- 
ment, long-lived stars and remnants, pair-instability supernovae, and black holes is tracked 
(see Table 1). The remnant neutron star or black hole mass is determined, when appli- 
cable, by Mr/bh = MAX(Mstar0.1^-'^-^^^®/^=*^% I.SSMq). This prescription roughly follows 
the results of Timmes et al. (1996) for population 111 stellar remnants. Metal enrichment 
is included using Z = 0.4(Mstar - 13.5MQ)/Mejccta (Woosley & Weaver 1982). Whenever a 
supernova occurs, for either population II or III, AfcjectalO^°erg Mq^ is deposited directly 
into the cell containing the star. 

We use a box size of 589 hr^ kpc, and the cosmology is set to h = 0.719, f2matter = 0.258, 
fl\ = 0.742, and Qh = 0.045. The power spectrum normalization is set to cxg = 1.0. The 
box size is smaller than recommended above, but it is intended to be a proof-of-concept 
simulation. The highest refinement at a given expansion factor a is kept at a physical 
resolution Axmin ~ 0.5pc, which corresponds to refinement level I — 17 at z ^ 13. For dark 
matter particles only, the highest level of refinement used in the cloud-in-cell approximation 
is / = 13 to avoid two-body relaxation effects (e.g., Knebe et al. 2000; Levine et al. 2008). 
The initial conditions were generated using the COSMICS package (Bertschinger 1995). The 
simulation box has a 256 cubed region (236 h~^ kpc) nested within a coarser 128 cubed grid. 
The dark matter particle mass in the higher resolution region is about 1000 Mq. 



5.2. Simulation Results 

Figure 3 shows maps of the average number density and temperature along the line of 
sight in a 360 pc cubed region at ^ ~ 13.3. The region contains a 5.3 x IO^Mq dark matter 
halo with about IO^Mq of gas within i?2oo ~ 420 pc. Two star clusters are present in this 
halo, with a combined mass ~ 9 x 10^ Mq. Most of this mass is in pop II stars, but about 
138Mo is from the pop III generation, either trapped in neutron stars or in stars that will 
become white dwarfs. Although pop III black holes were not formed in these clusters, but 
some clusters in the simulation have a total black hole mass > lOOOM©. More than half of 
the combined mass occurred in a burst of star formation that lasted about ~ 30 Myr. The 
effects of supernovae on the surrounding gas can be seen in the temperature and density 
plots, where bubbles have formed. The initial enrichment for the stars in Figure 3 was set by 
population III stars with a combined mass of ~ 19OM0. The final metallicity distribution is 
shown in the bottom-left panel, with a peak [Fe/H] ~ —1.3, consistent with Morgan & Lake 
(1989). Greif et al. (2007) found that a 200 M0 population III star, which did not occur in 
our simulation, could expel the gas from a halo, but that owes to the hard radiation from the 
precursor star and the extreme injection of 10^^ erg. The simulation demonstrates that the 
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formation of star clusters and their subsequent evolution can be captured from cosmological 
initial conditions. If these clusters survive to be called BGCs today, they represent a mode 
of GC formation that takes place in dark matter minihalos. This formation channel could 
operate until reionization limits cooling (Bromm & Clarke 2002). Cluster formation may 
then be delayed until massive disk galaxies form, which may be responsible for the red GC 
population. 

6. Conclusions 

During BGC formation, 3 x 10^ to 8 x IO^Mq of stars should have been produced. A large 
fraction, if not all, of the halo stars are remnants of destroyed BGCs and lower mass clusters 
that formed with BGCs. Using the conditional mass function, we found that the entire halo 
population could be produced by 2; ~ 13 with star formation efficiencies of ~ 10%. Based 
on constraints from reionization, radial biasing, and Press-Schechter theory, we argue that a 
simulation of a biased, 1 Mpc comoving region of the universe should capture the formation 
of BGCs by z ~ 10. If such a simulation does not capture the formation of these systems, 
it most likely reflects a problem in the subgrid model or cosmology. Finally, we presented 
a simulation that shows the formation of BGC-like clusters with masses ~ 5 x IO^Mq by 
z ~ 13. These clusters form in dark matter minihalos and can self-enrich to metallicities 
that are consistent with observed BGCs. 

This research was supported by an SNF Grant and the Zurich GTS. Simulations were 
run on the zbox machines, maintained by ITP UZH. We thank the anonymous referee for 
comments that improved this manuscript. We thank Ben Moore and Joachim Stadel for 
useful discussions and for making the Diemand et al. data available. 
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Mass Range Mq 


Outcome 


tsN (Myr) 


< 8 






8 < M < 25 


Mll, Mejecta 


10 


25 < M < 35 


Mbh, Mejecta 


10 


35 < M < 140 


Mbh 




140 < M < 260 


-Mejecta 


3 


260 < M 


Mbh 





Table 1: The mass in long-lived stars and remnants, Mll, black holes, Mbh, and ejecta mass 
Mejecta- The remnant mass is determined by Mr/bh = MAX(Mstar0.1^3-5^®/^=*-, I.35M0), 
which is based on the results of Timmes et al. (1996). The ejecta mass is the difference 
between the star's mass and its remnant. In the case of a pair-instability supernova (140 < 
M < 26OM0), the ejecta mass is the star's mass. The third column is the supernova delay 
time. 
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1. — The solid curve represents the MWG GCMF, as seen today, normahzed by 100 



clusters. The lines represent the ICMF mdN/dm oc m^^. The lowest line shows the os- 
culation mass mosc, the dashed line accounts for evaporation, and the dotted line includes 
evaporation+disk-shocking. 
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Fig. 2. — Conditional mass function for a halo of lO^^M© at z — 0. The curves represent the 
baryon fraction in collapsed objects for different dark matter halo minimum mass cutoffs. 
The thin, horizontal lines show the baryonic mass that is required to produce 5 x IO^A^q of 
halo stars and BGCs for a global star formation efficiency of 1% (top) and 10% (bottom). 
The thick line shows the fraction of mass that comes from v > i>2.5 (§4). The gray region 
shows the WMAP5 measurement for reionization. 
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[hi^H] 

Fig. 3. — Snapshot from a RAMSES simulation with nonequihbrium chemistry and the 
star formation prescription described in §5. The images correspond to z ~ 13.3. Top-left: 
average number density (cm~^) along the line-of-sight. Top-right: average density-weighted 
temperature over the mean molecular weight (T//i K). Bottom-left: star clusters in a dark 
matter halo ~ 5 x IO^Mq. Bottom-right: histogram of stellar metallicities for the stars 
shown in the bottom-left panel. 



